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We analyze a model of spinless fermions on a triangular lattice at half-filling interacting via 
strong nearest-neighbor repulsive interactions, V , using the variational Monte Carlo simulation 
technique. The existence of three-sublattice long-range order is confirmed by the finite-size 
scaling analysis of the charge structural factor at V c /t> 12. This ordered phase shows charac- 
teristics expected for a so called "pinball liquid" state, which has the spontaneous separation 
of fermionic degrees of freedom into coexisting Wigner crystal-like charge order (pin) and a 
metal (ball). The pins are fixed in order to maximize the kinetic energy gain of balls which 
move almost freely. The Fermi surface is reconstructed at V = V c and remains towards the 
strong coupling limit. These features reminiscent of the strong correlation together with the 
large value of V c /t distinguishes the pinball liquid from the conventional charge-density-wave. 

KEYWORDS: triangular lattice, honeycomb lattice, supersolid, organic solid, charge ordering, variational 
Monte Carlo method, geometrical frustration 



1. Introduction 

Geometrically frustrated lattices can provide par- 
ticular situations in which strong correlation destroys 
conventional classical orders and sometimes gives rise 
to another exotic states. Some simple examples are 
found in frustrated quantum spin systems, e.g. spin- 
gapped singlet states, 1_3 which are realized in realis- 
tic quasi low-dimensional materials such as CaV40g 4 
and SrCii2 (BC>3)2 - 5 ' 6 Another recent example is the ne- 
matic order which has been observed in spin systems 
on a square lattice with added ring exchange or frus- 
trated next-nearest-neighbor interactions. 7 ' 8 The quan- 
tum S = 1/2 spin systems are often described by the 
hard core bosonic model. The model is studied inten- 
sively and the supersolidity is realized under the frus- 
trated interactions on the triangular and the square lat- 
tices. 9 ~ 12 

The present paper deals with the spinless fermionic 
system on the triangular lattice. The basic nature of the 
ground state of this model is already studied by part of 
the authors; 13 ' 14 a partially charge ordered liquid called 
a "pinball liquid" is proposed based on the finite-size 
cluster calculation. 13 This state shows a three-sublattice 
structure, in which the carriers on one sublattice are 
localized as a Wigner-crystal (pin) and the remainders 
form a liquid (ball) itinerating along the other two sub- 
lattices. When we interpret the pin and ball as a solid 
and supcrfluid of the supersolid on the hard core bosonic 
model, we find that the pinball liquid has the same ori- 
gin to the supersolid, i.e., order by disorder mechanism. 13 
This suggests that although the statistics of bosons and 
fermions differ in two-dimension, the way they correlate 
is similar. 



However, comprehensive analysis of the fermionic sys- 
tem in two-dimension is far difficult from the case of 
bosons which is already known to have a long-range or- 
der in both the solid and the liquid part. 9 Therefore the 
details of the fermionic ground state nature still remain 
unclear and several questions arise, e.g. whether there 
is a long-range order of solid or not, what point differs 
from the conventional charge-density- wave (CDW), what 
kind of metallic state we have (Fermi liquid or non- Fermi 
liquid), or even whether the mctallicity remains stable 
in the bulk limit. To clarify these issues, we adopt the 
variational Mote Carlo (VMC) method to the half- filled 
ground state of the strongly interacting fermionic model 
on a triangular lattice. A charge correlation function in- 
dicates the realization of three-sublattice type long-range 
order, where two of the three sites arc equivalent. At the 
same time, the existence of the Fermi surface is suggested 
in the momentum distribution function. 

We shall mention that the topic is closely related to the 
experimental reports on organic solid, 6-ET2X, where 
a non-linear conductivity and many different kinds of 
short-range diffuse spots are reported in the low tem- 
perature metallic state. 15 ~ 17 This organic system is of- 
ten studied in terms of the extended Hubbard model 
(EHM) with interactions between charges populating on 
the same site (U) and on neighboring sites (V). 18 ~ 22 In 
such cases further difficulty arises in analyzing this model 
due to large degrees of freedom including spins. In fact, 
there remains some discrepancies in the phase diagram, 
e.g. the exact diagonalization studies show that the sys- 
tem is a "charge liquid" 18 without any long-range or- 
der in analogue to the spin liquid. In contrast, the VMC 
and the density matrix renormalization group (DMRG) 
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Fig. 1. The lattice structure and the transfer integrals t and the 
Coulomb repulsions V. 



study shows that the charge degrees of freedom has 
translationally broken structure in the metallic phase at 
V ~ V > 3£. 21 ' 22 We notice the fact that the strong cou- 
pling limit of the EHM without spin degrees of freedom 
can share the same classical nature as the present t-V 
model at half-filling we deal with. Therefore, our results 
can give a clue to understand the quantum fluctuation 
effect of charge degrees of freedom from that Ising limit 
of the EHM. Comparison of the results of the two models 
will be given in §4. 

2. Formulation 

The t-V model of spinlcss fcrmions is given explicitly 



as. 



n 



= E 



Tit + TLv 

(-tclcj+h-c. + Vmrij). (1) 



Here, cj (q) is the creation (annihilation) operator of 
fermions at site i and m = c\a. We allow a hop- 
ping between nearest-neighbors as well as next-nearest- 
neighboring sites along only one diagonal direction on 
a square lattice with a hopping amplitude t (t > 0) as 
shown in Fig. 1. The hopping integral t is taken as the 
unit of energy and V is the near-neighbor Coulomb re- 
pulsion energy. The sum on (i,j) runs over pairs of neigh- 
boring sites. Then, the dispersion is given by 

E(k) = -2t cos(k x ) - It cos(k y ) - 2i cos^ + k y ). (2) 

This model is equivalent to the equilateral triangular 
lattice. Here, we remark for the half-filled case in or- 
der to clarify the features of the pinball liquid state in 
the strongest geometrical frustration. The similar analy- 
sis for other filling, focusing on Na2Co04, was reported 
in ref. 23. 

To investigate ground state with strong interactions, 
we compute the variational energy, E = (^\H\^) / (^\^) , 
by using the VMC method. Since the Hamiltonian con- 
sists only of the near-neighbor coupling term, it is nec- 
essary to add a Jastrow 24 term to the wave function to 
lower the variational energy: we study the intersite cor- 
relation effect by using Jastrow-type trial wave function 
|\&t) as defined by 



I*t) =Ph 



(3) 



where \<p) is an one-body wave function determined by 
diagonalizing the mean-field Hamiltonian. Ph is a near- 



neighbor projection operator given by 

(ij) 



(4) 



where ft is a variational parameter which controls the 
weight of the wave function; when the fermion exists on 
both i and j sites, the wave function is multiplied by 
ft, while it remains unchanged for the other cases. In the 
limit of ft. = 1, the wave function corresponds to \<p) itself, 
while at ft — > the occupation on nearest-neighbor site 
is completely suppressed. 

In the present paper, we adopt the L x L (L = 6 ~ 18) 
lattices with periodic boundary conditions in both di- 
rections, and optimize the energy expectation values for 
different V with the total Monte Carlo step greater than 
3 x 10 7 . The system does not satisfy the closed shell con- 
dition for half-filling. However, we confirmed that there is 
no qualitative difference in the result of closed shell with 
periodic-antiperiodic boundary condition. We choose two 
different kinds of one-body wave functions, |</>) = |</>fs) 
and |</>cdw), which are the Fermi sea and the CDW solu- 
tion of the mean-field calculation, respectively. The for- 
mer has a translational symmetry and gives the simplest 
and unbiased results when projected. We also examine 
the possibility of the breaking of the translational sym- 
metry by adopting the latter CDW wave function. The 
comparison of these two cases gives us useful informa- 
tion. 

3. Results 

3.1 Projected Fermi sea 

First, we assume the projected Fermi sea, P/i |0fs)> as 
a variational wave function for the normal phase. The In- 
dependence of the optimized variational parameter ft opt 
is shown in Fig. 2. The ft op t rapidly decreases from the 
non- interacting value, ft op t = 1, as V increases. The func- 
tional form changes at V ~ 12, which suggests a change 
of the character of the state. The expectation values of 
kinetic part Ti t of eq. (1) is also shown in the inset as a 
function V. It shows a nearly 1^-independent insulating 
behavior at V > 12, however, even in this region there is 
a finite kinetic energy gain. 

Next, to investigate the instability towards charge or- 
der, we measure the charge correlation function by the 
optimal wave function, which is defined as 



C(k) 






(n 3 ))c 



i(vi—rj )-k 



(5) 



where N s (= L x L) is the number of sites. In inset 
of Fig. 3(a), two peaks are observed in C(k) at k = 
±(2-7r/3, 27r/3) for V = 12. This wave number corre- 
sponds to the charge modulation patterns based on a 
three-sublattice structure. As shown in Fig. 3(a), the 
peak in C(k) grows as V increases and starts to satu- 
rates above the inflexion point, V > 10. 

To see whether the state has a long-range order, 
we perform a finite-size scaling analysis and examine 
whether the charge structure factor C(Q) diverges in the 
thermodynamic limit. If the long-range order exists, the 
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Fig. 2. The optimized near-neighbor correlation factor /i op t as a 
function V by using the Jastrow type wave function based on the 
Fermi sea. The system is 12 X 12. The inset shows the expectation 
value -Eopt of Tit part as a function of V. 



charge structure factor is expected to behave as 

C(Q) = N S A 2 + O (v^\C) , (6) 

where we expect C(Q)/N S to follow a linear function of 
l/\/N s with a nonzero intercept which gives the bulk 
order parameter A. We plot C(Q)/N S versus l/\/W s in 
Fig. 3(b). The nonzero extrapolation value appears when 
V > 12, which implies the existence of a long-range order 
based on the three-sublattice structure above the phase 
transition point, V c ~ 12. 25 

We shall emphasize that in this scaling analysis the 
wave function has translational symmetry, so that V c 
gives a reliable value in the bulk limit. Actually, V c has 
the same order of magnitudes as the onset value of the 
supersolid order in the hard core bosonic t-V model by 
the quantum Monte Carlo study, 9 t/V c ~ 0.12, and is 
consistent. If we consider instead the projected CDW 
wave function with three- fold periodicity (i.e., transla- 
tional symmetry is broken by hand), the value of Vcdw 
shall be several times smaller as we show in the next 
section. 

3.2 Projected charge- density-wave 

In the presence of the long-range order, the transla- 
tional symmetry is broken in the bulk limit. However, 
in the present finite-size cluster, the projected Fermi sea 
wave function has the translational symmetry, and the 
orders appear as the correlation which lasts towards the 
bulk limit. There shall instead be a symmetry- broken 
wave function that could give better description to the 
bulk charge ordering. As such a candidate, we employ an- 
other projected wave function, P^I^cdw); where |</>cdw) 
is obtained as a CDW solution of a mean-field Hamilto- 



Hmf = y^(-fcfcj + h.c. 



(i,j) 



Y^ (Ai cos (Q ■ n) + A 2 sin (Q ■r i ))i 



(7) 
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Fig. 3. (a) Intensity of the peak of the charge correlation func- 
tion C(fe) at k = Q = ±(27r/3,27r/3) as a function of V. The 
trial wave function, Ph\<f>~Fs), is adopted on the N s = 12 X 12 
lattice. The inset shows the fc-dependence of C(fc) at V = 12. 
(b) System-size dependence of charge structure factors C(Q) for 
V = 12 and V = 24. The dotted lines are least-squares fits to 
the data. 



where Ai and A2 are CDW order parameters and Q = 
±(27r/3,27r/3) is a wave vector which characterizes the 
state. 

Charge configurations of the mean-field solutions ex- 
hibit three different classes, A-B-B, A-A-B and A-B-C 
type, as shown schematically in Figs. 4, where A, B, and 
C denote the charge-rich sublattice, the charge-poor sub- 
lattice, and the sublattice with the intermediate charge 
density (nearly charge-neutral), respectively. The param- 
eters (Ai > 0, A 2 = 0) and (Ai < A 2 = 0) correspond 
to A-A-B and A-B-B type of structure, respectively. The 
parameters (Ai = and A 2 ^ 0) realize an A-B-C 
type structure. In the VMC calculation, we take h and 
Ai(or A 2 ) as variational parameters, and find the opti- 
mized solution with finite Ai > at V > Vcdw ~ 3. 
This A-A-B type solution has lower energy than the 
A-B-C type one. The energy difference is estimated as 
(£abc - E AAB )/N S = 0.014 at V = 10 on a N = 12 x 12 
cluster. This result is different from the exact diagonal- 
ization study which shows A-B-C type of charge correla- 
tions, 13 which is possibly attributed to the smaller clus- 
ter size in that previous study. In the following, only the 
charge configuration of the A-A-B type is considered. 

Figure 5 shows the V-dependcnce of h opt of the pro- 
jected A-A-B type wave function. Following the emer- 
gence of finite Ai > at V > Vcdw ~ 3, h opt shows sig- 
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Fig. 4. Representative charge ordered patterns studied in this pa- 
per. Charge density of (a) A-B-B type, (b) A-A-B type, and (c) 
A-B-C type. Black and white circles represent charge-rich A and 
charge-poor B sites, respectively. Gray circles denote intermedi- 
ate C sites. 
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5. Optimized /i op t as a function of V based on the CDW trial 
wave function, P/, | <j>cuw ) ■ '^' ne rGsu lt of fhe Fermi sea function, 
Ph\<t>¥s)i (Fig. 2) is plotted here together for comparison. Inset 
shows the expectation value of charge density along the x-axis 
on the N s = 12 X 12 lattices at V = 24 by the Phl^cow)- 



nificant deviation and the energy is lowered compared to 
those of the projected Fermi sea. The expectation value 
of the charge density, pi = (qcj), is given in the inset of 
Fig. 5. The amplitude of pins on B-sitc is almost 0, which 
means that the holes are strongly localized as a pin. Such 
a hole-pinned configuration is expected to have the same 
energy with the electron-pinned configuration (i.e. A-B-B 
type with pa = 1) at half- filling if one assumes a particle- 
hole (p-h) symmetry. In the present fermionic system we 
always have the hole-pinned state, since the triangular 
lattice has p-h asymmetry. If we take the negative hop- 
ping amplitude t < 0, an electron-pinned pinball liquid 
is expected instead. 

The rigorous pinball liquid state is created just off the 
filling factors, p = 1/3 and 2/3, by adding a single elec- 
tron or a hole, respectively, at t/V — > 0. 13 Combined 
with the present results, we may conclude that the pin- 
ball liquid state is stable over the entire range of the 
charge density 1/3 < p < 2/3. This phase is further 
classified into a hole-pinned (p ~ 2/3) and an electron- 
pinned (p ~ 1/3) ones, and the transition between them 
shall exists at p > 1/2. 
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Fig. 6. Momentum distribution functions for various V on the 
N s = 12 X 12 system; Panel (a) and (b) show the ones based on 
PjJ^pg) and Phl'f'cvw^' res P ec ti ve ly; which are plotted within 
the first Brillouin zone of the normal state, i.e. along the arrows 
in panel (c). The inset of panel (b) shows h(k) at V = 24, which is 
the replotted n(k) within the folded Brillouin zone for the three- 
sublattice, i.e. along the arrows in panel (d). Panels (c) and (d) 
show the Fermi surfaces of an non-interacting fermions on the 
triangular lattice at half-filling and on the honeycomb lattice at 
1/4-filling, respectively. In (d) the Brillouin zone folded by 1/3 
from the one in (c), which is shown in broken line. The dots on 
panel (d) denote the representative data points of the momentum 
distribution function in the inset of panel (b). 



3.3 Fermi surface reconstruction 

In order to understand the origin of the breaking of 
translational symmetry in more detail, we analyze the 
momentum distribution function at the Fermi surface, 



n(fe) 



-y> 



c\ Cj )e 



i(ri—rj)-k 



(8) 



Figures 6(a) and 6(b) show the ^-dependence of n(k) 
measured with the optimized parameters along the path 
(0, 0)-(0, 7r)-(7r, 7t)-(0, 0) in the Brillouin zone for the trial 
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wave functions, Ph\4>Fs) and Ph\<t>cDw) > respectively. For 
comparison, we show the Fermi surface for the non- 
interacting case in Fig. 6(c). 

In Fig. 6(a) we find a large jump- like structures in n(k) 
for Ph\(pFs) a t V = 0, 3, 7. The positions of the jumps are 
consistent with the Fermi surface of the non-interacting 
system. When V = 24, the jump is blured and instead 
a structure at around (w, tt) appears. Thus, it is likely 
that the long-range order at V = 24 is related to the 
reconstruction of the Fermi surface. 

Actually, the reconstruction distinctively takes place 



when we adopt the wave function P; 



'h\ < PcBWi 



at V = 24; 

new discontinuities are found both along (0, 7t)-(7T, tt) and 
along (tt, 7r)-(0, 0) in Fig. 6(b). It is detected more clearly 
by assuming the three-fold periodicity and considering a 
momentum distribution function defined on the reduced 
Brillouin zone, 



ra(fc) = ^n(fc + G), 



(9) 



G 



with G being the reciprocal lattice vector, as shown in 
the inset of Fig. 6(b) along the path (0, 0)-(0, 2tt/3) and 
(7r/3,7r/3)-(0, 0). In Fig. 6(d), we also plot the represen- 
tative values of n(k) within the folded Brillouin zone. 

Let us focus only on A-sitcs in Fig. 4(b). Namely, 
when the hole-pinned B-sites are depleted from the tri- 
angular lattice, remaining carriers propagate on a hon- 
eycomb lattice formed by A-sites. In Fig. 6(d), we also 
plot the Fermi surface of the non-interacting system on 
the honeycomb lattice at 1/4-filling. Note that the hon- 
eycomb bands have a double structure within the first 
Brillouin zone, and at 1/4-filling (which corresponds to 
iVs/6-balls on 2iV s /3-sites) there exist a half- filled band 
and an empty band so that the large ellipsoidal Fermi 
surface covers half the Brillouin zone. The present result 
shows that the location of the discontinuity oih(k) possi- 
bly falls on the honeycomb Fermi surface, which suggests 
that this Fermi surface describes well the metallic state 
of the pinball liquid at V > V c . 

Therefore, it seems that the large Fermi surface at V = 
(Fig. 6(c)) is discontinuously transformed to another 
large Fermi surface at V > V c . At the transition, the 
number of carriers decreases by N s /3 which corresponds 
to the number of fixed pins. 

The ellipsoidal Fermi surface in the non-interacting 
case suggests that the charge susceptibility xo only takes 
a broad peak structure 26 disadvantageous of the nesting 
instability. Therefore, we expect to have large fluctua- 
tions of charges at V — > V c over a certain range of wave 
numbers around the broad peak. Even in such case, the 
mean-field solutions may allow for a phase transition into 
the CDW state which is usually considered as an ar- 
tifact of underestimating the charge fluctuation. Since 
Ph\(f>GDw) assumes the periodicity by hand in a mean- 
field manner, it also overestimates the three-fold struc- 
ture although the true long-range order is still absent. 

In fact, our results show that the sign of reconstruction 
is already found at V — 3, while we still find a distinct 
jump near (0, tt) in n(k) in both V = 3, 7 similar to the 
non-interacting case, as depicted in the main panel of 
Fig. 6(b). Therefore it is reasonable to consider that the 



Fermi surface at V < V c is close to the non-interacting 
one in Fig. 6(c). The precursor of the pinball liquid struc- 
ture appears below the true long-range order transition 
due to the short range correlation, while the true recon- 
struction of the Fermi surface may take place only at 
V >V C . The fact that V c is four times as large as Vcdw 
also suggests that these two cannot merge and that the 
Fermi surface reconstruction at V c is not attributed to 
the nesting instability. 

At V = 24, the correlation effect is as strong as to have 
^opt ~ 0.2 -c 1, while h(k) (inset of Fig. 6(b)) behaves 
almost like a step function of the free fcrmions. This may 
indicate that once the pins arc fixed the balls move al- 
most freely in a one-body- like manner in this phase which 
extends towards the strong coupling limit. 

3-4 Energy gain 

In order to characterize the phase above the transi- 
tion, we finally analyze the charge condensation energy, 
AE = E(P h \cf> FS }) - E(P h \^^)), which is the en- 
ergy gain in breaking the translational symmetry The 
kinetic and the interaction part of the condensation en- 
ergy are separately plotted as a function of V in Fig. 7. 
For comparison, we also show the condensation energy 
of the mean-field wave function in the inset, AE = 
E(Ph\<f>Fs)) — E(\(J)mf))- This purely mean-field state is 
obtained by setting h = 1, i.e.,|</>MF) = I'/'cdwH' 1 = !)• 
The CDW transition takes place at different points, 
Vcdw ~ 3 and 2 for the projected and purely mean- 
field wave functions, respectively, reflecting the degree 
of correlation included. Therefore, we shall again stress 
that the CDW phase is artifactual and here we discuss 
it just for references. 

In the CDW phase, 3 < V < 9, in order to break the 
translational symmetry one shall cost the kinetic energy 
while gaining the interaction energy, i.e., AEy > and 
AE t < 0. This mechanism is essentially the same as the 
case of mean- field wave function (see the inset of Fig. 7), 
which is a weak coupling picture. This suggests that the 
obtained projected CDW state is a mean-field like state, 
and by including the correlation effect in a more proper 
manner in the presence of translational symmetry, the 
kinetic energy gain shall be enhanced. 

On the other hand, at V > 10, the kinetic energy 
gain is maximized while there is a cost of the interac- 
tion energy, i.e., AEy < and AE t > 0, when adopt- 
ing the projected CDW wave function. It means that 
although both wave functions have long-range order, the 
explicit breaking of translational symmetry (projected 
CDW wave function) is energetically favored due to large 
E t ] in the projected CDW state, fixing the pins stabilizes 
the metalic path (balls are less scattered), i.e., kinetic en- 
ergy gain is maximized. In the projected Fermi sea the 
charges can fluctuate at most to avoid the interaction 
(the strong correlation effect). Such fluctuation is sup- 
pressed in the projected CDW state to some extent by 
the periodicity, which leads to the interaction energy loss. 
Thus, maximizing the kinetic energy gain plays a key 
role in stabilizing the long-range ordered pinball liquid 
at V > V c . 

Similar characteristic behavior of AEy and AE t at 
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Fig. 7. The kinetic and interaction parts of the optimized charge 
condensation energy, AE t = E t (P h \(pFs)) - E t(Ph\<t>cQw)) and 
AE V = E v (P h \c/) FS )) - Bv(-Pfc|^cDw))> respectively as a func- 
tion of V. The data is from 12 X 12 lattices. The inset shows the 
result for mean-field calculation. 



around the phase transition point is found by Yokoyama 
et al. in the superconducting phase transition and order- 
disorder transition in the Hubbard model, 27 which might 
be a universal features found in the strongly correlated 
systems. 

4. Summary and discussion 

In the present paper we studied the ground state of the 
two-dimensional spinless fermion t-V model on a triangu- 
lar lattice at half-filling by the VMC method. We found 
the strong indication for the long-range order (sponta- 
neous breaking of the translational symmetry) of the pin- 
ball liquid state at V > V c ~ 12; partially localized holes 
(not the electrons) on one of the three-sublattices form a 
long-range order, while the other two sublattices equiv- 
alently contribute to the metallic conduction, i.e., an A- 
A-B structure which breaks the particle-hole symmetry. 
It is stabilized by suppressing the fluctuation of pins, al- 
lowing for a stable (almost free) metallic path to gain 
kinetic energy, which is consistent with the strong cou- 
pling picture. 13 By the finite-size scaling analysis based 
on the wave functions with the translational symmetry, 
the results give the first evidence of the existence of the 
long-range order and of the metallicity stabilized towards 
the bulk limit at V > V c . 

Regarding the phase at Vcdw < V < V c , we do not 
have a conclusive answer. The three-fold periodic long- 
range order is absent in the bulk limit. Therefore, the 
CDW wave function obtained variationally by project- 
ing the mean-field solution does not give the best de- 
scription of this phase since the translational symmetry 
is broken by hand. The search for another variational 
wave function without symmetry breaking is beyond the 
present scheme. Meanwhile our results suggest a short- 
range order of fluctuating three-sublattice type of struc- 
ture. Therefore, we may safely conclude that the short- 
range order of charges grows with increasing V due to 
strong correlation effect and at V = V c , a long-range or- 
der appears. The reconstruction from a large Fermi sur- 
face to another large Fermi surface shall takes place at 



V = V c , where the number of carriers decreases by the 
number of crystallized pins. 

Let us briefly discuss the comparison with the re- 
lated theoretical studies of the electronic system on the 
anisotropic triangular lattice. The metallic charge order 
is found in the EHM at quarter-filling, which is to be 
compared to the present t-V model at half-filling; the 
exact diagonalization study mentions the idea of "charge 
liquid" without long-range order which completely over- 
looked the geometrically favored three-sublattice con- 
figuration due to small cluster size. 19 We stress that 
the frustration induced metallic state they discussed has 
long-range order and is to be discussed in analogy with 
the "supcrsolid" 9-11 and not the "spin liquid". Actually, 
it is figured out in the VMC study on EHM based on 
the mean-field wave function as the A-B-B type of order 
at Vc ^ 4, while the DMRG gives A-A-B type of order 
a t V c > 3. The reason for the different types of configu- 
ration in the same EHM model is unclear. However, the 
DMRG gives the reasonable hole density, n^ ~ 1.5 and 
Ub ~ 0, which is well regarded as the pinball liquid but 
not as the CDW (see the inset of Fig. 5). 

As for the value of their V c (without scaling analysis) 
four times smaller than our estimation we speculate two 
possible reasons; one is that in the EHM the pins con- 
sists of the doubly occupied holes, 22 and at (7 < oo, 
it may help to stabilize the long-range order. The mag- 
netism which is incompatible with the pinball liquid is 
suppressed by this double occupancy. Another possibil- 
ity is that our Vcdw ~ 3 based on the mean-field CDW 
wave function corresponds to their V c , which is actually 
close to the mean-field value. Therefore, in the former 
VMC study, the same treatment with the present study 
might give quantitatively the similar result. The latter 
DMRG might also overestimate the orders due to Fricdel 
oscillations under the open boundary condition, and the 
true long-range order may shift to larger values of V. 

Regarding the short-range orders in the weak coupling 
phase, the DMRG study is consistent with our interpre- 
taion. 22 The competition of three-sublattice and stripe 
type of charge ordering is discussed in more detail in the 
RPA study by including the electron-phonon(e-h) cou- 
pling terms to the EHM. 26 The phase transition they 
discuss corresponds to our Vcdw an d we also consider 
that it may remain as short range order if one excludes 
the effect of e-h coupling. From all these reports, we shall 
conclude that in both t-V model and EHM, a short-range 
order exists in the weak coupling region which transforms 
into the long-range ordered pinball liquid at strong cou- 
pling although there remains some quantitative discrep- 
ancies. The robustness of the metallicity is guaranteed 
even in the presence of spin degrees of freedom and for 
systems with different types of non-interacting density of 
states. 

The detailed comparison of the above mentioned theo- 
retical findings and the experiments on organic 0-ET2A 
systems 15-17 are given elsewhere, 21 ' 26 and we just pro- 
vide in this paper the essential picture of the charge 
degrees of freedom. Finally, we shall mention that not 
only in organic solids but in the transition metal oxide, 
AgNiC-2 28 there exists a charge ordered metallic state 
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which is possibly interpreted as a pinball liquid state. 
Although there is an orbital degrees of freedom in this 
material, the essential picture should be the same. We 
expect that more materials will be found that show such 
particular features of the strongly correlated metallic sys- 
tem on the similar geometrically frustrated lattices. 
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